Effect of Biomolecular Conformation on Docking Simulation: A Case Study on a Potent HIV-1 Protease Inhibitor.

Human immunodeficiency virus infection/acquired immunodeficiency syndrome (HIV/AIDS) is a disease pertained to the human immune system. Given its crucial role in viral replication, HIV-1 protease (HIV-1 PR) is a prime therapeutic target in AIDS therapy. In this regard, the dynamic aspects of ligand-enzyme interactions may indicate an important role of conformational variability in HIV-1 PR inhibitor/drug design. In the present contribution, the effect of HIV-1 PR flexibility (within multiple crystallographic structures of HIV-1 PR) on binding to the Amprenavir was elucidated via an ensemble docking approach. Molecular docking studies were performed via advanced AutoDock4.2 software. Ensemble docking of Amprenavir into the active site of various conformations of HIV-1 PR predicted different interaction modes/energies. Analysis of binding factors in terms of docking false negatives/positives revealed a determinant role of enzyme conformational variation in prediction of optimum induced fit (PDB ID: 1HPV). The outcomes of this study demonstrated that conformation of receptor may significantly affect the accuracy of docking/binding results in structure-based rational design of anti HIV-1 PR agents. Furthermore; some strategies to re-score the docking results in HIV-1 PR targeted docking studies were proposed.


Introduction
Acquired immunodeficiency syndrome (AIDS) is a disease related to the human immune system (1). Human immunodeficiency virus (HIV) has been identified as the etiological agent of AIDS (2). Cells of the immune system, called T-cells or CD4 cells that are responsible for fighting against infections and other physiological disturbances are attacked and destroyed by HIV. One of the essential HIV enzymes, whose activity is necessary for viral replication, is HIV-1 protease (HIV-1 PR) (3). In fact, production of mature and infectious viral particles is depended on the proteolytic activity of the HIV-1 PR and for this reason; this enzyme was recognized as a major therapeutic target in AIDS therapy and structures of HIV-1 proteases (17,18). Amprenavir ( Figure 1) is a potent and selective HIV-1 PR inhibitor with sub-nanomolar HIV-1 PR inhibition activity (k i =0.6 nM) (19,20) and hence was selected as a model in our studies.
HIV-1 PR inhibitors are believed to inactivate the HIV-1 protease leading to the immature, non-infectious viral particles (6). Most of the developed HIV-1 protease inhibitors are peptidomimetic molecules (7). The main drawback of peptidomimetic compounds is their low oral bioavailability arising from high molecular weight and poor solubility (8). Due to this limitation, many researchers have focused on nonpeptidic HIV-1 PR inhibitors (3,9,10).
Amprenavir, Atazanavir, Darunavir, Indinavir, Fosamprenavir, Lopinavir, Nelfinavir, Ritonavir, Saquinavir and Tipranavir are typical anti-AIDS drugs that have been approved by the United States Food and Drug Administration (US FDA) as HIV-1 PR inhibitors. These drugs are currently used in combination therapy with reverse transcriptase inhibitors (11,12). Although several successful drugs have been developed against AIDS, current status shows a rapid emergence of drug resistance to most of the HIV-1 PR inhibitors (11). In this regard, recent research aimed at proposing new anti-protease agents with minimum side effects and being able to delay the appearance of resistance (13,14).
In continuation to our interest in structure based modeling of bioactive molecules (15,16) and to further elucidate the important role of target conformation in molecular docking results, we decided to explore the significance of HIV-1 PR flexibility through ensemble docking of Amprenavir into the multiple crystallographic into related carbon atoms of the receptor and Kollman charges were also assigned. For docked ligands, non-polar hydrogens were merged; Gasteiger charges assigned and torsions degrees of freedom were also allocated by ADT program. 100 independent genetic algorithm (GA) runs were considered. 2.5×10 7 maximum number of evaluations was used for Lamarckian GA method. All other docking parameters were set at their default values. A grid of 60×60×60 points in x, y, and z direction was built centered on the center of mass of the catalytic site of HIV-1 PR crystallographic structures. Cluster analysis was performed on the docked results using a root mean square (RMS) tolerance of 2 A˚.
Schematic 2D representations of the ligandreceptor interactions were all generated using LIGPLOT (23).

Docking validation
A performance of a docking simulation method was checked via its ability in reproducing a binding mode for a cocrystallographic (cognate) ligand (24). For this purpose, the structure of a cognate ligand (Amprenavir) was retrieved and re-docked into the active site of HIV-1 PR structures. Root mean square deviations (RMSD) of the Cartesian coordinates of the re-docked ligand atoms proved the validation of docking method for further modeling studies (Table 1) (25). As it is obvious from the summarized data, all the crystallographic files under study represented adaptable predictability level (26) within 100 independent genetic algorithm (GA) runs and 2.5×10 7 maximum number of evaluations for Lamarckian GA method. It should be noted that those structures exhibiting RMSD values over 3 may also pass the filter when considering their number of active torsions (27).

Ensemble docking of Amprenavir
We aimed to evaluate the Amprenavir / HIV-1 PR interaction considering ligand induced enzyme conformation. Our dataset included one apo and fifty holo HIV-1 PR structures. These structures were subjected to ensemble docking procedure. Crystallographic structure of the Amprenavir/HIV-1 PR complex was deposited in the PDB website (1HPV) (28) and as mentioned before, this crystallographic structure was considered as the reference point in our docking simulations.
The RMSD of the backbone carbon atoms (Cα) in the selected PDB structures ranged 0.22-0.85 and 0.24-0.93 Å in chains A and B of HIV-1 PR, respectively (with regard to the PDB code: IHPV; Figure 2). Different RMSD values indicated the conformational changes of HIV-1 PR upon binding to the various inhibitors.

Amprenavir/HIV-1 PR interactions
Lipophilic contacts ( Figure 4) and H-bond interactions (Table 3) in docked Amprenavirprotease complexes were monitored. According to the 2D Ligplot diagrams, thirty-two residues of the HIV-1 PR were found to make lipophilic contacts with Amprenavir within fifty-one enzyme conformational structures. In the case of hydrogen bond interactions, a total of nineteen amino acids interacted Amprenavir within 51 conformations of the enzyme. Data are summarized in Table 3 while numbers refer to the H-bond distances.

Validation of virtual binding affinities
To further validate the AutoDock binding affinities, two co-crystallographic HIV-1 PR/ inhibitor datasets with available biological activities at PDB bind (29) (Figure 5) or Binding MOAD (30) ( Figure 6) databases were selected for a regression analysis. AutoDock binding affinities were all obtained from the self-docking step.

Effective factors in binding conformation
We were interested in finding the factors that might be determinant in induced conformation of HIV-1 PR/Amprenavir complex (PDB ID: 1HPV). For these purpose; a binding system comprised of three major constituents (ligand, enzyme and their interaction) was taken into consideration. Such a system may be defined by several descriptors that are related to the system constituents ( Figure 7).
To account for the conformational deviation of HIV-1 PR from its apo structure (native conformation), a pair wise structure alignment     Table 4.

Analysis of binding results via docking false negative/positives
Ensemble docking approach may be interpreted in terms of predicted false negative (FN)/false positive (FP) results. High rate of false negatives/positives is a common issue in docking procedure leading to low "hit rates". Due to this rationale, we decided to evaluate the docking results (Table 4) via FP and FN results.
For the sake of clarity, estimated descriptors (binding factors) for Amprenavir/HIV-1 PR co-crystallographic complex (IHPV) were considered as reference points in our analysis. In this manner, two distinct regions may be considered for each binding factor; a distance between reference level and optimum level including FPs and a distance between reference level and non-optimum level including FNs ( Figure 8).

Ensemble docking approach
Docking is a popular virtual structure-based method that is used in the design of biologically interesting molecules (31). It enables the          prediction of stereoelectronic complementary fit of a potential bioactive ligand with its biomolecular target. In this regard; availability of crystallographic data on HIV-1 PR (Brookhaven protein databank website: http://www.rcsb. org) facilitated the performance of structure based drug discovery projects aiming at HIV-1 PR as a biomolecular target for AIDS disease.

Code of docked
The HIV-1 PR consists of two identical 99 amino acid monomers representing a homodimer with C2 symmetry. Each subunit includes one of the two conserved triads (Asp-Thr-Gly) containing the catalytically active aspartate residues; Asp 25 and Asp 25′ (32). It has been well known that upon binding of different HIV-1 PR inhibitors, significant conformational changes might be expected for the enzyme (33, 34). Indeed, dynamic aspects of binding in the interaction of HIV-1 PR inhibitors with HIV-1 PR active site are crucially important for the design of novel enzyme inhibitors. Due to the computational cost in designating numerous degrees of freedom, incorporation of meaningful protein flexibility during a docking procedure is a difficult task although several efforts have been performed (35).
One of the alternative approaches for the flexible-receptor docking is the cross-docking of a typical ligand into the multiple crystallographic structures of the receptor (protein ensemble structures) (21). Holo crystallographic structures of targets provide appropriate models that represent real ligand induced conformations upon binding to the various chemical scaffolds (different inhibitors). In the case of biological targets lacking sufficient crystallographic holo structures, conformational ensemble may be generated virtually. However the advantage of the latter approach would be the possibility of generating more protein conformations but at the same time, a major drawback remains; the produced protein conformations may not be indicative of real structures.
A simple flow chart representing the ensemble docking procedure might be depicted as below (Figure 9). It should be noted that ligand binding ensembles (resulted from ensemble docking) may be subsequently exploited as valuable input data for quantitative structure binding relationship studies.
Results of ensemble docking showed that Amprenavir interacted with HIV-1 PR active site via different binding modes. None of the docked ligands showed completely identical binding poses in the active site of the HIV-1 PR and the best scored conformation might not be supported with highest binding energy (refer to Table 2).

Amprenavir/HIV-1 PR Interactions
The frequency of occurrence for a specific chemical interaction in multi-conformational ligand-enzyme assemblies may indicate the significance of such interaction in ligandenzyme complex. Regarding the binding data, some principles might be driven: -Docked Amprenavir showed different hydrophobic and H-bond binding patterns in multiple conformational ensembles of HIV-1 PR active site.  (Table 3). Analysis of binding maps showed that hydroxyl group of Amprenavir contributed to the H-bond(s) with Asp25(A) and Asp25(B) while sulfonamide oxygen atoms may be involved in H-bond interactions with Ile50(A) and Ile50(B) residues of HIV-1 PR. 2D schematic representation of binding interactions between Amprenavir and HIV-1 PR structure (PDB ID: 4DJO) is depicted in Figure 11.

Regression analysis of docking results versus biological data
Our regression analysis showed that docking outputs could be used for the elucidation of HIV-1 PR inhibitory activities (PDB bind database) with a relatively good predictability level (R 2 =0.703; Figure 5). Our results exhibited a lower regression coefficient for Binding MOAD activities (R 2 =0.443; Figure 6).

Analysis of binding determinants
We decided to rank the probable determinant factors of the ligand induced enzyme conformation. In our opinion, the results of such study might assist in re-scoring the docking results within a screened dataset. Results of pair wise alignment study with apo conformation of the enzyme (3IXO) showed that co-crystallographic HIV-1 PR/Amprenavir complex (PDB ID: 1HPV) was associated with minimum geometrical deviation of enzyme from its apo structure (RMSD=0.38 Ǻ, Table 4). This observation confirmed the literature evidence that in binding to the inhibitors, a majority of enzymes might be necessarily redecorated via an optimum geometrical path (36). However, the most geometric deviation of the enzyme could be observed for the HIV-1 PR conformation designated by PDB code 2XL2 (RMSD=0.87 Ǻ, Table 4). For further consideration, 2D schematic representation of pair wise structural alignments between chains A of apo HIV-1 PR (3IXO) and holo HIV-1 PRs (IHPV and 1XL2) were depicted in Figure 12. Analysis of residues showed that maximum distortion of HIV-1 PR conformation in 1XL2 structure occurred within a loop containing Gly48, Gly49, Ile50, Gly51, Gly52 and Phe53 residues (red highlighted in Figure 12b). Estimated descriptors for various Amprenavir/HIV-1 PR systems (Table 4) were normalized (0-100%) to elucidate their probable significance in achieving the optimum target conformation upon binding to Amprenavir (PDB ID: 1HPV). For this purpose, each descriptor was designated by two numerical values indicating the optimum and non-optimum levels. In the   i.e., N-benzyl-2-(2,6-dimethylphenoxy)-N-[((3R,4S)-4 {[isobutyl(phenylsulfonyl) amino] methyl} pyrrolidin-3-yl)methyl] acetamide, the most distorted residues are highlighted by red circle in 1XL2. case of ∆E instability of docked ligand conformation and RMSD of enzyme from reference structure, generally accepted optimum values were 0 kcal/ mol and 0 Ǻ, respectively, hence these values were taken as optimum levels. It should be noted that no commonly accepted thresholds for optimum scores of AutoDock binding affinity, number of lipophilic interacted residues and number of H-bond interactions in a typical enzyme/inhibitor system could be rationalized. Due to this restriction, optimum levels of these factors were considered as the best achieved scores within the docked Amprenavir/HIV-1 PR systems. Similarly, the lowest numerical levels of all descriptors were taken as the worst achieved scores within the docked Amprenavir/HIV-1 PR systems (Table 5). The probable significances of five descriptors in the achieved induced fit of Amprenavir/HIV-1 PR complex (1HPV) were reported as significance percentages. Significance values of the descriptors were all estimated within the optimum and non-optimum levels ( Table 5).
Data mining showed that induced fit of Amprenavir/HIV-1 PR complex might be significantly determined by lipophilic contacts followed by deviation of enzyme from its native conformation, H-bond patterns, estimated free binding energy and deviation of ligand from its optimum conformation (designated by ∆E instability ), respectively. Of course we believe that such priority order have been achieved within the selected dataset in this study and more extended explorations through larger enzyme/inhibitor datasets would be less biased to the size of dataset.

FP and FN
Analysis of binding factors exhibited that none of the Amprenavir conformations could be recognized as FP points on the basis of factor d (conformational variation of enzyme from apo structure). This observation is very important and emphasizes on the determinant role of enzyme conformational variation in prediction of ligand induced binding poses. Further investigations via chemically diverse inhibitors may be possibly the subject of future investigations in this field.
There was an opposite case for factor e (conformational variation of ligand from optimum structure); forty-two FP points could be predicted. Such a result may be translated into the uncertainty of factor e in prediction of HIV-1 PR targeted docking results and confirmed our previous results that inhibitors might not necessarily interact with the enzyme active site via their minimum energy conformation (18,19). This was also in agreement with our above analysis on binding factors i.e., significance percentage of 21.4% was estimated for factor e ( Table 5).
AutoDock binding affinities (factor a) and number of H-bond interactions (factor c) produced relatively balanced results (Figure 8). However analysis of docking results on the basis of hydrogen binding exhibited twenty-two non-FP/FN points. It should be noted that most of the H-bond patterns showed reasonable agreement with the binding pattern of Amprenavir in its crystallographic file (IHPV).
Most of the Amprenavir conformations were predicted as FNs on the basis of lipophilic interactions (seven FPs and thirty-seven FNs) but less non-FP/FNs were resulted (5 points).
The outcomes of this study revealed that a major problem in docking based virtual screening is the proper selection of an enzyme conformation. Following this rationale and on the basis of results taken form ensemble docking approach, different scenarios may be considered: 1) Docking validation (self-docking) protocols may be performed with less trouble due to the presence of induced target structure.
2) Our ensemble docking approach on HIV-1 PR system demonstrated that varied binding results might be expected upon docking of a specific inhibitor (Amprenavir) into the multiple conformations of the enzyme. To alleviate the problem, a simple docking approach within an enzyme including a similar cognate (cocrystallographic) ligand (similar holo structure) followed by an efficient scoring function is proposed.
3) In the case of holo enzyme structures bearing non-similar cognate ligands, an ensemble docking approach may be run by the cross-docking of a co-crystallographic enzyme inhibitor into the multiple enzyme structures (holo dataset). Subsequent analysis of probable induced fit determinants (section 3.1) may be done within the results of ensemble docking approach. Ranked induced fit determinants could be used in post-scoring of the ensemble docking results. 4) When no holo structure is available, an ensemble docking approach may be run through apo structures of the enzyme.

Conclusion
Computer aided molecular design (CAMD) has spurred a renewed interest to deal with the growing body of information from genomic and proteomic efforts. In this regard, molecular docking is an attractive branch of CAMD that allows drug designers to simulate binding mode and predict binding affinity of different ligandreceptor complexes. In the present study, ensemble docking approach was successfully applied for modeling of anti-AIDS agent Amprenavir in the active site of HIV-1 PR. The outcomes of this study showed that success of a typical HIV-1 PR targeted docking strategy in rational drug design might be strictly depended on a selection of docked enzyme conformation. Further results showed that in selection of a desirable HIV-1 PR target for docking of amprenavir like ligands, lipophilic contacts are very important while the effect of ligand departure from its optimum conformation is less important. Pertaining to this, the multiple-receptors docking approach might be a suitable strategy to find a relatively optimum conformation of the enzyme to run the docking simulation of a query class of inhibitors. It is apparently known that our analysis method might be biased due to the restricted dataset of crystallographic files, but retrieved protein conformations (PDB database) can be regarded as valuable sources of such studies since they represent real induced enzyme conformations upon binding to the assayed inhibitors. Moreover; the results of ensemble docking approach may be complementary to molecular dynamics simulations and hence assist in finding optimum dynamic paths.